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We study the transitions from band insulator to metal to Mott insulator in the ionic Hubbard 
model on a two dimensional square lattice using determinant Quantum Monte Carlo. Evaluation 
of the temperature dependence of the conductivity demonstrates that the metallic region extends 
for a finite range of interaction values. The Mott phase at strong coupling is accompanied by 
antiferromagnetic (AF) order. Inclusion of these intersite correlations changes the phase diagram 
qualitatively compared to dynamical mean field theory. 
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Introduction 



Interaction effects in tight-binding models such as the 
Hubbard Hamiltonian have been widely studied, and un- 
derstood, for their ability to drive transitions to magneti- 
cally ordered states and insulating behavior. Also exten- 
sively studied, but less well understood, is the converse 
phenomenon, namely the possibility that correlations can 
cause metallic behavior. For this latter problem, at- 
tention has focused on an Anderson insulating starting 
point, in which the electrons are localized by random- 
ness. The primitive picture is that, especially away from 
commensurate fillings where correlations tend to localize 
particles in Mott insulating (MI) states, electron-electron 
repulsion can spread out the wave function and cause 
delocalization. Despite considerable effort, whether this 
effect actually occurs for fermionic particles in two di- 
mensions, that is, whether metallic phases exist for dis- 
ordered, interacting 2d electron systems, is not settled 1 . 
In an attempt to gain leverage on the problem, the phases 
of disordered and interacting bosonic systems have now 
been studied. Models such as the boson Hubbard Hamil- 
tonian arc more amenable to exact numerical studies and 
much is now known about the insulating, glassy, and su- 
perfluid phases, and transitions between them£. 

A somewhat simpler context in which to study the pos- 
sibility of interaction driven insulator to metal transitions 
is to begin with a band insulating (BI) state, in which 
the insulating behavior is caused by a periodic external 
potential as opposed to a random one&i Recently, this 
issue has been addressed within dynamical mean field 
theory (DMFT) and a number of interesting conclusions 
emerged 5 . However, because DMFT treats only a single 
site (retaining, however, all the dynamical fluctuations 
of the self-energy which is ignored in conventional, static 
mean field theory) it is important to undertake comple- 
mentary work which is able to retain intersite fluctua- 
tions. 

In this paper, we investigate such Bl-metal transitions 
with determinant Quantum Monte Carlo (DQMC). We 



study the "ionic Hubbard model" : 

H = - t {c] a ci a + c\ a c 3a ) + U "/T«U 
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where cj (c^-) are the usual fermion creation (destruc- 
tion) operators for spin a on site I, and ni a — c[ a ci a is the 
number operator, t, fi and U are the electron hopping, 
chemical potential, and on-site interaction strength, re- 
spectively. The kinetic energy sum is over near neighbor 
sites (Ij) on a two dimensional square lattice. A(— 1) ; 
is a staggered site energy. In the noninteracting limit, 
U = 0, the e ffect of A is to produce a dispersion relation, 
E(k) = ±y/e(k) 2 + A 2 with e(k) = -2t[cosk x + cosk y ], 
which is gapped at half-filling. A considerable amount is 
known concerning this model in one dimension*, but the 
existence of an interaction driven metallic phase at half- 
filling is still unresolved even in d = 1. Metal-insulator 
transitions in a related system with randomly located site 
energies with a bimodal distribution have also been stud- 
ied within DMFTifi. 

In this letter we will use DQMC to study the role of 
interactions in driving a Bl-metal transition in the model 
described by Eq. JIJ. 



Computational Methods 

DQMC 9 provides an exact numerical approach to 
study tight binding Hamiltonians like the Hubbard 
model. The partition function Z is first expressed as 
a path integral by discretizing the inverse temperature 
/3. The on-site interaction is then replaced by a sum over 
a discrete Hubbard-Stratonovich fieldiS. The resulting 
quadratic form in the fermion operators can be integrated 
analytically leaving an expression for Z in terms of a sum 
over all configurations of the Hubbard-Stratonovich field 
with a summand (Boltzmann weight) which is the prod- 
uct of the determinants of two matrices (one for spin up 
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and one for spin down). The sum is sampled stochasti- 
cally using the Metropolis algorithm. The results capture 
correlations in the Hubbard Hamiltonian exactly since 
the systematic 'Trotter errors' associated with the dis- 
cretization of the inverse temperature can easily be ex- 
trapolated to zero. Results must also be extrapolated to 
the thermodynamic limit, as we shall discuss^. 

Equal time operators such as the density and energy 
are measured by accumulating appropriate elements, and 
products of elements, of the inverse of the matrix whose 
determinant gives the Boltzmann weight. We will show 
results for the spin structure factor, 

S(k) = Xy k l ( ( n J+'T ~ n i+u)( n iT - n ji) ) > 
I 

which probes magnetic order. For the conductivity, <7d c , 
We employ an approximate procedure^ which allows 
Ode to be computed from the wavevector q- and imag- 
inary time r-dependent current-current correlation func- 
tion h xx {<\, t) without the necessity of performing an an- 
alytic continuation^, 

2 

CT dc = — A xx (q = 0,r = /3/2) . 

IT 

Here (3 = 1/T, A xx (q,r) = (j x (q, r)j x (-q, 0)), and 
jx (q, t) the q, r-dependent current in the x-direction, is 
the Fourier transform of, 
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This approach has been extensively tested and used for 
the superconducting-insulator transition in the attractive 
Hubbard modeliS, as well as for metal-insulator transi- 
tions in the repulsive modeli4. 
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FIG. 1: The transitions, at half-filling, from a band insulator 
to metal to MI with increasing U axe shown for periodic po- 
tential strength A = 0.5. At (7=0 the conductivity <Tdc goes 
to zero as T is lowered. However, for U = It, 2t the system 
is metallic. Mott insulating behavior sets in for U = At. The 
lattice size is 6 x 6. When A = 1.0, the band gap increases 
and U = It is no longer sufficiently large to screen the one 
body potential and drive the system metallic. 




Results 

We begin by showing the temperature dependence of 
the conductivity <7dc for increasing values of the inter- 
action strength for A = 0.5. In Fig. ^ we see that the 
insulating behavior at U — 0, signaled by da^c/dT > 
at low T, is changed to metallic dadc/dT < at low T 
when U — 1. A further increase of the correlations to 
U = 2 weakens the metallic behavior, which is finally 
destroyed completely in a transition to a MI at U = 4. 
When the band gap is larger (A = 1), the screening of the 
one-body potential is not sufficiently strong for U = 1 to 
cause metallic behavior, as is shown by the correspond- 
ing data set in Fig. ^ Unless otherwise mentioned, the 
lattice size used in the simulations is N = 6 x 6 and the 
filling is p = 1.0 (half-filling). 

In the single site (t = 0) limit, the ionic Hubbard model 
is a band insulator for U < 2A and a MI for U > 2A. 
That is, at weak coupling and half-filling, the sites with 
lower energy —A are doubly occupied and those with 



FIG. 2: The conductivity a^ c at half-filling for A = 0.5 is 
shown as a function of U for three different low tempera- 
tures, j3 — 10, 12, 16. The band-insulator to metal transition 
is signaled by the crossing of the curves at U c i ~ OAt. At 
U C 2 ~ 2At the three curves cross again, indicating the MI 
transition. 



higher energy +A are empty, with a gap to further addi- 
tion of particle set by 2A — U. At strong coupling, both 
types of sites are singly occupied, with a 'Mott' gap to 
further addition of particles set by U — 2 A. At the single 
special value U = 2A correlations close the gaps^. Fig- 
ure |2 which presents results for Ode for A = 0.5, shows 
that when t is nonzero, this single metallic point is ex- 
panded to a finite range of U values. Interestingly, how- 
ever, the largest conductivity remains near U = 2A = 1 
as one might expect from the t — analysis. The BI to 
metal transition occurs at U c \ « OAt, where the change 
in the order of the three curves indicates a transition 
from (Tdc decreasing as (3 increases to tide increasing as (3 
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FIG. 3: The AF structure factor is shown at half-filling as a 
function of U for A = 0.5 and = 10, 12, 16. 
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FIG. 4: The conductivity <7dc is shown as a function of tem- 
perature at half-filling. When the periodic potential, and 
hence the non-interacting band gap, is absent (A = 0.0) 
the square lattice Hubbard model is insulating for all U, 
due to nesting of the Fermi surface. We re-display data for 
A = 0.5, U — 1 from Fig.Qto emphasize the contrast between 
the metallic behavior there and the insulating behavior for all 
U when A = 0. 



increases. The metal to MI transition is at U c i w 2At, 
where <7d c once again decreases as (3 increases. 

The use of DQMC to study the ionic Hubbard model 
allows us to examine the behavior of intersite correla- 
tions, such as the spin-spin correlations and their Fourier 
tranform 5(k). Fig. [21 shows results for the AF structure 
factor S(ir, tt) as a function of U for (3 — 10, 12, 16. Com- 
paring with Fig. |21 we see that the band insulating and 
metallic phases are paramagnetic, but that the transition 
to MI behavior is accompanied by the onset of AF order. 

One way in which the inclusion of such intersite corre- 
lations changes the physics in a fundamental way is when 
the periodic potential is absent, that is, at A = 0. In 
DMFT in the paramagnetic phase, the Hubbard model 
is a metal at weak couplingi^i^. However it is well known 
that the d = 2 half-filled square lattice Hubbard model, 
Eq. ^ is an AF insulator at all U, even weak coupling. 




FIG. 5: The conductivity at half-filling is shown for different 
lattice sizes for (7 = 1, close to the point where the system 
is most metallic for periodic potential A = 0.5. (See Fig. [2]) 
Although <7d c decreases with increasing lattice sizes, the sig- 
nature of metallic behavior (dadc/dT < 0) is unchanged. 



Figure 01 presents our results for the conductivity which 
confirm this. At all U values shown, Cdc ultimately de- 
creases as T is lowered. Indeed, we have verified that 
the value of T where a^ c has its maximum correlates 
well with the temperature T* at which AF correlations 
begin to rise rapidly. This temperature, like the Neel 
temperature in the d = 3 Hubbard model, is a non- 
monotonic function of U, fa lling to small values both at 
weak (T* oc texp(—ay/t/U) and at strong (T» oc t 2 /U) 
coupling. To our knowledge, this is the first time the 
insulating nature of the square lattice Hubbard model 
at weak coupling has been shown from Quantum Monte 
Carlo studies of Ode- It is interesting to note that while 
all the A = curves share a common low temperature 
slope dadc/dT > 0, a distinction between the origins of 
insulating behavior in the AF and Mott is clearly evi- 
dent. At small U , (Tdc attains a much larger value before 
turning over as T is lowered than in the strong coupling 
Mott regime. 



While DQMC allows us to look at intersite corre- 
lations and concomitant phenomena like antiferromag- 
netism, the method employs lattices of finite size, unlike 
DMFT which directly probes the thermodynamic limit. 
Thus, it is important to verify that the metallic phase 
we observe persists on larger lattices. In Fig. we show 
results for Ode as a function of temperature in the metal- 
lic phase for lattices up to 12 x 12. The rise in a^c with 
decreasing T is seen to occur for all the lattices stud- 
ied. We comment that it is not surprising that we find 
the lattice size has a rather substantial influence on the 
conductivity for these parameters, since it is known that 
such finite size effects are larger at weak coupling. 
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FIG. 6: The phase diagram of the ionic Hubbard model. Sym- 
bols are the result of our QMC simulations. The dashed line 
is the strong coupling (t=0) phase boundary between band- 
insulator and Mott insulator. 



Conclusions 

We have presented determinant Quantum Monte Carlo 
studies of the two-dimensional ionic Hubbard Hamilto- 
nian which demonstrate that interactions can drive a 
band insulator metallic. This work complements DMFT 
studies by including intersite AF correlations which qual- 
itatively alter the ground state phase diagram. 

We have focused most of our results on A = 0.5. 
However, we have also performed simulations sweeping 
U at A = 0.25 and A = 1.00. The emerging phase 



diagram is shown in Fig. 6. There are several key dif- 
ferences with that obtained with DMFT—. First, as we 
have emphasized, the behavior along the A = axis is 
significantly altered. Contrary to DMFT, the inclusion 
of intersite magnetic fluctuations yields an AF insulat- 
ing phase for all U. Like the DMFT treatment, however, 
we find a metallic phase intervening between band and 
Mott insulators for nonzero A. This phase is centered 
roughly around the strong coupling boundary (dashed 
line in Fig. 6). 

The sign problem prevents us from performing simu- 
lations much beyond A/t ks 1. However, we expect that 
the sign problem will become better in the limit of large 
A, where we have very widely separated bands. Related 
studies of the boson-Hubbard model in a "superlattice" 
potential, which exhibit a band-insulator to superfluid 
transition^, show the appearance of insulating phases 
at half-integer fillings. These 'charge-transfer' insulators 
occur as a result of Mott splitting of the widely separated 
bands ' . We plan to explore this possibility in the fcrmion 
case in future work. 

We acknowledge support from NSF DMR 0312261, and 
useful input from S.A. Clock. 

Note Added: A recent preprint 17 reports simulation re- 
sults for the ionic Hubbard model using cluster DMFT. 
Like our approach, this method incorporates intersite 
correlations, and as found here, obtains a Mott phase 
along the entire A = axis of the phase diagram. A key 
difference is that the cluster DMFT approach suggests 
the intermediate phase between Mott and band insula- 
tors is a correlated 'bond ordered' insulator. 
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